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A TECHNIQUE OF TREATING NEGATIVE WEIGHTS IN WENO SCHEMES 

JING SHlt, CHANGQING HU*, AND CHI- WANG SHU§ 


Abstract. High order accurate weighted essentially non-oscillatory (WENO) schemes have recently been 
developed for finite difference and finite volume methods both in structured and m unstructured meshes 
A key idea in WENO scheme is a linear combination of lower order fluxes or reconstructions to obtain a 
higher order approximation The combination coefficients, also called linear weights, are determined by 
local geometry of the mesh and order of accuracy and maj become negative WENO procedures cannot 
be applied directly to obtain a stable scheme if negative linear weights are present Previous strategy for 
handling this difficulty is by either regrouping of stencils or reducing the order of accuracy to get rid of the 
negative linear weights In this paper we present a simple and effective technique for handling negative linear 
weights without a need to get rid of them Test cases are shown to illustrate the stability and accuracy of 
this approach 

Key words, weighted essentially non-oscillatory, negative weights, stability, high order accuracy, shock 
calculation 

Subject classification. Applied and Numerical Mathematics 

1. Introduction. High order accurate weighted essentially non-oscillatory (WENO) schemes have re- 
cently been developed to solve a hyperbolic conservation law 

u t + V f(u) = 0 (11) 

The first WENO scheme was constructed in [18] for a third order finite \olume version in one space dimension 
In [10], third and fifth order finite difference WENO schemes in multi space dimensions are constructed, with 
a general framework for the design of the smoothness indicators and nonlinear weights Later, second, third 
and fourth order finite volume WENO schemes for 2D general triangulation have been developed in [4] and 
[8] Very high order finite difference WENO schemes (for orders between 7 and 13) have been developed in 
[1] Central WENO schemes have been developed in [12], [13] and [14] 

WENO schemes are designed based on the successful ENO schemes in [7, 23, 24] Both ENO and 
WENO use the idea of adaptive stencils in the reconstruction procedure based on the local smoothness 
of the numerical solution to automatically achieve high order accuracy and non-oscillatory property near 
discontinuities ENO uses just one (optimal in some sense) out of many candidate stencils when doing the 
reconstruction, while WENO uses a convex combination of all the candidate stencils, each being assigned 
a nonlinear weight which depends on the local smoothness of the numerical solution based on that stencil 
WENO improves upon ENO in robustness, better smoothness of fluxes, better steady state convergence, 
better provable convergence properties, and more efficiency For a detailed review of ENO and WENO 
schemes, we refer to the lecture notes [21, 22] 
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WENO schemes have already been widely used in applications Some of the examples include dynamical 
response of a stellar atmosphere to pressure perturbations [3], shock vortex interactions and other gas 
dynamics problems [5], [6], incompressible flow problems [26], Hamilton- Jacobi equations [9], magneto- 
hydrodynamics [11], underwater blast- wave focusing [15], the composite schemes and shallow water equations 
[16], [17], real gas computations [19], wave propagation using Fey’s method of transport [20], etc 

A key idea in WENO schemes is a linear combination of lower order fluxes or reconstructions to obtain 
a higher order approximation The combination coefficients, also called linear weights, are determined by 
local geometry of the mesh and order of accuracy and may become negative WENO procedures cannot 
be applied directly to obtain a stable scheme if negative linear weights are present Previous strategy for 
handling this difficulty is by either regrouping of stencils (e g m [8]) or reducing the order of accuracy (e g 
in [12]) to get rid of the negative linear weights In this paper we present a simple and effectne technique 
for handling negative linear weights without a need to get rid of them Test cases will be shown to illustrate 
the stability and accuracy of this approach 

We first summarize the general WENO reconstruction procedure, consisting of the following steps We 
assume we have a given cell A (which could be an interval in ID, a rectangle in a 2D tensor product mesh, 
or a triangle in a 2D unstructured mesh) and a fixed point x G within or on one edge of the cell 

1 We identify several stencils S 3 , j = 1, ,q, such that A belongs to each stencil We denote In 

Q 

T = U Sj ^ ie larger stencil which contains all the cells from the q stencils 
j=i 

2 We have a (relatively) lower order reconstruction or interpolation function (usually a polynomial), 
denoted by Pj(x), associated with each of the stencils S } , for j = 1, ,q We also have a (relatively) 
higher order reconstruction or interpolation function (again usually a polynomial), denoted by Q(x), 
associated with the larger stencil T 

3 We find the combination coefficients, also called linear weights, denoted by 71, , j q , such that 

Q(x G ) = ^2 'Y j Pj{x g ) (12) 

j=i 

for all possible given data in the stencils These linear weights depend on the mesh geometry, the 
point x G , and the specific reconstruction or interpolation requirements, but not on the given solution 
data in the stencils 

4 We compute the smoothness indicator, denoted by /?,, for each stencil Sj, which measures how 
smooth the function Pj(x) is in the target cell A The smaller this smoothness indicator 0j, the 
smoother the function Pj(x) is in the target cell In all of the current WENO schemes we are using 
the following smoothness indicator 

&= £ f \^ a[ ~\D a Pj {x)Ydx (13) 

i<M<t Ja 

for j = 1, ,q, where k is the degree of the polynomial Pj (x) and |A| is the area of the cell A in 2D 

This factor is different for ID or 3D the purpose of it is to bring the smoothness indicator invariant 
under spatial scaling 

5 We compute the nonlinear weights based on the smoothness indicators 


w, 


CJ, 


7 


(e + 0j) 2 


( 14 ) 
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where 7 , are the linear weights determined in step 3 above, and e is a small number to avoid the 
denominator to become 0 We are using e = 10 -6 in all the computations in this paper The final 
WENO approximation or reconstruction is then given by 

R ( x °) = YL^jPjix 0 ) (15) 

j=i 

We remark that all the coefficients in the above steps which depend on the mesh but not on the data of 
the numerical solution, should be computed and stored at the beginning of the code after the generation of 
the mesh but before the time evolution starts 

We now use a simple example to illustrate the steps outlined above We assume we are given a uniform 
mesh I, = {x,-i/ 2 ,x l+ i/ 2 ) and cell averages of a function u(x) in these cells, denoted by u , We would 
like to find a fifth order WENO reconstruction to the point value u{x l+ \/ 2 ), based on a stencil of five cells 
, /,+i , 4 + 2 }, with the target cell containing the point x l+ i/ 2 chosen as A — I, 

In step 1 above we could have the following three stencils 

£l = $2 = {4-1, In 4+i}, S 3 = {I l ,I t+1 ,I t+2 }, 

which make up a larger stencil 


T — {1,-2, I t -i , I,, I t +i, It+ 2 } 


In step 2 above we would have three polynomials p 3 (x) of degree at most two, with their cell averages 
agreeing with that of the function u in the three cells in each stencil S 3 The higher order function Q{x) is a 
polynomial of degree at most four, with its cell averages agreeing with that of the function u in the five cells 
in the larger stencil T The three lower order approximations to u(x l+ 1 / 2 ), associated with Pj(x), in terms 
of the given cell averages of u, are given by 

, x 1 - 7_ 11_ 

£ 1 ( 27 + 1 / 2 ) = 2 > Ul ~ 2 ~ 6 U ' _1 + 

, . 1 _ 5 _ 1 _ ,, 

P2{X l+ i/ 2 ) = ~Q U t~ 1 + + 3 U i+l> l 1 G ) 

, X 1- 5 _ 1_ 

P3\X 2) — + gW«+l qUi+2 

Each of them is a third order approximation to u(x t+] / 2 ) The higher order approximation to u(x l+ j/o), 
associated with Q(x), is given by 


, 1 13 47 9 1 

Q ( x i + i / 2 ) - 35^-2 - ~ 

which is a fifth order approximation to u(x, +1 /o) 

I 11 step 3 above we would have 


1 3 

7i - 10 > 72 - 5 , 



It can be readily verified, using (1 6 ) and (1 7), that 


( 17 ) 


< 3 ( 27 + 1 / 2 ) = 71 Pi ( 27 + 1 / 2 ) +72 P2(x, + 1 / 2 ) + 73 £ 3 ( 27 + 1 / 2 ) 
for all possible given data u 3 , j = i — 2,i — 1, 1 , t + 1 , 1 + 2 
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Fig 1 1 Reconstructions to u(x t+ 1 / 2 ) Solid lines exact function, symbols numerical approximations Left fifth order 
WENO Right fifth order traditional 


In step 4 above we could easily work out from (1 3) the three smoothness indicators given by 

13 _ , l _ 2 

A = — (u,_ 2 - 211,-1 + u,)- + - {u ,- 2 - 4u,_i + 3 u,)~ , 

13 2 l 2 

fh = {ui - 1 - 2 u, + u,+i) + - (u t - 1 - u l+ 1)‘ , 

13 _ 2 1 _ - 2 

/?3 = — (u t — 2u, + i + u l+ •>)' + - (3u, — 4u, + i + u,+ 2 )“ 

Iz 4 

We notice in particular that the linear weights 71 , 70 , 73 in step 3 above are all positive In such cases, the 
WENO reconstruction procedure outlined above and the scheme based on it work very well In Fig 11 we 
plot the approximation to u{x) for a discontinuous function u(x) = 2x for x < 0 and u{x) = —20 otherwise, 
by the fifth order WENO reconstruction on the left and by the fifth order traditional reconstruction (17) 
on the right, with a mesh x t = (2 — 0 4965)Ax with Ax = 0 02 We can clearly see that WENO avoids the 
over and undershoots near the discontinuity 

We now look at another simple example where some of the linear weights in step 3 above would become 
negative We have exactly the same setting as above except now we seek the reconstruction not at the cell 
boundary but at the cell center x, This is needed by the central schemes with staggered grids [12] Thus, 
step 1 would stay the same as above, step 2 w r ould produce 

1 _ 1 _ 23. 

P\\ x i) — 24 U, ~ 2 i2 U,_1 24 Ul ' 

, , 1 _ 13 _ 1 „ n . 

P2 ^ Xl ^ = _ 24 U, -i + 12 U! “ 24 fi *+i» ( 1S ) 

. , 23 _ 1 1 

P3( x i) = Y\ u i + 12 U,+1 _ 24 u,+2 

Each of them is a third order reconstruction to u{x t ) The higher order reconstruction to u(x t ), associated 
with Q(x), is given by 

, 3 29 1067 _ 29 _ 3 

Q{X,) ~ 640 U,_2 “ 480 Ut_1 + 960 U 1 ~ 48Q Ui+1 + 640 U,+2 ’ 


4 


( 19 ) 




Fig 1 2 Reconstructions to u(x t ) Solid lines exact function, symbols numerical approximations Left fifth order 
WENO Right fifth order traditional 


which is a fifth order reconstruction to u(x t ) Step 3 would produce the following weights 

9 49 9 

7i-_ 80’ 72 ~ 40’ 73 ~ _ 80 

Notice that two of them are negative The smoothness indicators in step 4 will remain the same This 
time, the WENO approximation, shown at the left of Fig 1 2, is less satisfactory (in fact, even worse than 
a traditional fifth order reconstruction show on the right), because of the negative linear weights 

We remark that negative linear weights do not appear in finite difference WENO schemes m any spatial 
dimensions for conservation laws for any order of accuracy [10], [1], and they do not appear in one dimensional 
as well as some multi-dimensional finite volume WENO schemes for conservation laws Unfortunately, the\ 
do appear in some other cases, such as the central WENO schemes using staggered meshes we have seen 
above, high order finite volume schemes for two dimensions described in [8] and in this paper, and finite 
difference WENO approximations for second derivatives 

While on approximation alone the appearance of negative linear weights might be annoying but perhaps 
not fatal (Fig 1 2), in solving a PDE the result might be more serious As an example, in Fig 1 3 we show 
the results of using a fourth order finite volume WENO scheme [8] on a non-umform triangular mesh shown 
at the left, which has negative linear weights, for solving the two dimensional Burgers equation 



(1 10 ) 


in the domain [—2,2] x [—2,2] with an initial condition uo(x,y) = 0 3 + 0 7sin ( j(x + y)) and periodic 
boundary conditions We can see that serious oscillation appears in the numerical solution once the shock 
has developed The oscillation eventually leads to instability and blowing up of the numerical solution for 
this example 

The mam purpose of this paper is to develop a simple and effective technique for handling negative linear 
weights without a need to get rid of them Test cases will be shown to illustrate the stability and accuracy 
of this approach 
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Fig. 1.3. 2D Burgers’ equation. Left: non-uniform triangular mesh used in the computation. Right: fourth order WEND 
result at t — 0.473, cfl=0.2, without any special treatment for the negative linear weights. 


2. A splitting technique. We now introduce a splitting technique to treat the negative weights. It 
is very simple, involves little additional cost, yet is quite effective. The WENO procedure outlined in the 
previous section is only modified in step 5 in the following way: 

5’ If mini'll,..., y q ) > 0 proceed as before. Otherwise, we split the linear weights into two parts: 
positive and negative. Define 

it = \ (7/ + 0 )7i|) , 7f = 7 f - 7 i = 1 , (I (2.1) 

where we take 8 = 3 all the numerical tests. We then scale them by 

a± = E?f ; 7* = 7? /a*, . * = 1, (2-2) 

j=l 

We now have two split polynomials 

Q ± (x g ) = jrif Pj (x G ) (2-3) 

i=i 

which satisfy 

Q(x G ) =- <j + Q + (x g ) - a ~Q~( x g ). (2.4) 

We can then define the nonlinear weights (1.4) for the positive and negative groups 7 * separately, 
denoted by ojf, based on the same smoothness indicator 0j. We will then define the WENO 
approximation J r ? ± (ar G? ) separately by (1.5), using and form the final WENO approximation by 

R{x g ) — <j + R + (x°) — cr~ Rr(x G ). 

We remark that the key idea of this decomposition is to make sure that every stencil has a significant 
representation in both the positive and the negative weight groups. Within each group, the WENO idea of 
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FlG. 2.1. WENO approximations with the splitting treatment, for negative linear weights. Left: approximation to u(x-;). 
Bight: Burgers equation , solution at t = 5/7 r 2 , cfl=0.2. 


redistributing the weights subject to a fixed sum according to the smoothness of the approximation is still 
followed as before. 

For the simple example of fifth order WENO reconstruction to u(xj), the split linear weights correspond- 
ing to (2.1) are, before the scaling, 

9 9 . 49 49 9 __ 9 

7l ~ 80’ 7l ~ 40’ 77 ~ 20’ 72 “ 40’ 7:i “ 80’ 7:5 ~ 40' 

We notice that, as the most expensive part of the WENO procedure, namely the computation of the 
smoothness indicators (1.3), has not changed, the extra cost of this positive/negative weight splitting is very 
small. 

However this simple and inexpensive change makes a big difference to the computations. In Fig. 2.1 we 
show the result of the two previous unsatisfactory cases, the fifth order WENO reconstruction to u(xi) in 
Fig. 1.2 left, and the approximation to the Burgers equation in Fig. 1.3 right, now using WENO schemes 
with this splitting treatment. We can see clearly that the results are now as good as one would get from 
WENO schemes having only positive linear weights. 

It is easy to prove that the splitting maintains the accuracy of the approximation in smooth regions. 
We will demonstrate this fact in the following sections. We will also demonstrate the effectiveness of this 
simple splitting technique through a few selected numerical examples in the next sections. The main WENO 
schemes we will consider are fifth order finite volume WENO schemes on Cartesian meshes, and the third 
and fourth order finite volume WENO schemes on triangular meshes. In both cases negative linear weights 
appear regularly. 

The calculations are performed on SUN Ultra workstations and also on the IBM SP parallel computer 
at TCASCV of Brown University. The parallel efficiency of the method is excellent (more than 90%). 

3. 2D finite volume WENO schemes on Cartesian meshes. 

3.1. The schemes. We describe two different ways to construct fifth order finite volume WENO 
schemes on Cartesian meshes. Comparing with finite difference WENO methods [10], finite volume meth- 
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ods have the advantage of an applicability of using arbitrary non-uniform meshes, at the price of increased 
computational cost [2] 

We define the cell 


h,j = [a\_ $,*,+ $] x [y J _i,y J+ i] (3 1) 

for i — 1, ,m,j = 1, ,n, Where I ltJ needs not be uniform or smooth varying 

The three-point Gaussian quadrature rule is used at each cell edge when evaluating the numerical flux 
in order to maintain fifth order accuracy Let (x G ,y G ) denote one of the Gaussian quadrature points at the 
cell boundary of I tJ given by T = {x = x ( _i,i/ J _i < y < y J+ 1 } There are two ways to perform a WENO 
reconstruction at the point (x G ,y G ) 

Genuine 2D The first WENO reconstruction is genuine 2D finite volume We can see that there are 
totally nine stencils S St t ( s,t = —1,0,1) Each stencil S s ,t contains 3x3 cells centered around I l+StJ+t 
On each stencil we can construct a Q 2 polynomial (tensor product of second order polynomials in x and y) 
satisfying the cell average condition (l e its cell average in each cell inside the stencil equals to the given 

l 

value) Let T = (J which contains 5x5 cells centered around I tJ On T we can construct a Q x 

s,t =- 1 

polynomial satisfying the cell average condition The WENO reconstruction is then performed according to 
the steps outlined in sections 1 and 2 

We would like to make the following remarks 

1 By using a Lagrange interpolation basis, we can easily find the unique linear weights 

2 Even for a uniform mesh, a negative linear weight appears for the middle Gaussian point (x G , y G ) = 
( x,_i, y, ) Such appearance of negative linear weights has also been observed in the central WENO 
schemes [12], see the example in sections 1 and 2 before 

3 By Taylor expansions, we can prove that the smoothness indicators yield a uniform fifth order 
accuracy in smooth regions See [10] for the method of proof 

Dimension by Dimension The second WENO reconstruction exploits the tensor product nature of the 
interpolation w’e use This WENO procedure is performed on a dimension by dimension fashion The WENO 
schemes applied in [5], [6] belong to this class Consider the point ( x G ,y G ) as above First we perform a 
one dimensional WENO reconstruction in the y direction, in order to get the one dimensional cell averages 
(in the x direction) w(»,y G ) Then w r e perform another one dimensional WENO reconstruction to w in the 
x direction, to obtain the final reconstructed point value at (x G ,y G ) 

We would like to make the following remarks 

1 For a scalar equation, the underlying linear reconstructions of the above two versions are equivalent 
For nonlinear WENO reconstructions they are not equivalent Both of them are fifth order accurate 
but the actual errors on the same mesh may be different, see Table 3 1 below 

2 For systems of conservation laws such as the Euler equations of gas dynamics, both versions of the 
WENO reconstruction should be performed in local characteristic fields 

3 The dimension by dimension version of the WENO reconstruction is less expensive and requires 
smaller memory than the genuine two dimensional version The CPU time saving is about a factor 
of 4 for the Euler equations in our implementation The computed results are mostly similar from 
both versions 

In the following, we will give numerical examples computed by the above WENO schemes Splitting 
technique has been used in all the computations when negative linear weights appear W T e will show the 
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Table 3 l 
2D vortex evolution 



genuine FV 

dim-by-dim 

N 

Ax 

L°° error 

order 

L°° error 

order 

20 

G 71E-1 

4 38E-2 


5 26E-2 


40 

3 77E-1 

3 10E-3 

4 59 

5 66E-3 

3 86 

80 

2 01E-1 

1 20E-4 

5 15 

3 96E-4 

4 22 

ICO 

1 00E-1 

4 39E-G 

4 7G 

7 96E-G 

5G2 

320 

5 00E-2 

1 88E-7 

4 53 

2 90E-7 

4 77 


results for both smooth and discontinuous problems 

3.2. 2D vortex evolution. First, i\e check the accuracy of the WENO schemes constructed above 
The two dimensional vortex evolution problem [21], [8] is used as a test problem 

We solve the Euler equations for compressible flow m 2D 

U t + f(U) x + g(U) y = 0, (3 2) 

where 

U = ( p,pu,pv,E ) T , 
f(U) = ( pu,pu 2 + p,puv,u(E+p)) T , 
g(U) = {pv, puv, pv 2 + p, v(E + p)) T 

Here p is the density, ( u , v ) is the velocity, E is the total energy, p is the pressure, related to the total energy 
by E — + \p{u 2 + v 2 ) with 7 = 14 

The setup of the problem is the mean flow is p = 1, p = 1, {u,v) = (1,1) and the computational 
domain is [0, 10] x [0, 10] We add, to the mean flow, an lsentropic vortex (perturbations in (u, v) and the 
temperature T = 2, no perturbation in the entropy S — 

(Su,Sv) = ~e 05 ^- r2 \-y,x), ST = 55 = 0 

where ( x,y ) = (a: — 5, y — 5), r 2 = x 2 + y 2 , and the vortex strength e = 5 

We use non-uniform meshes which are obtained by an independent random shifting of each point from a 
uniform mesh in each direction within 30% of the mesh sizes The solution is computed up to t = 2 Table 3 1 
shows the L°° errors of p We can see that both the genuine Wo dimensional finite volume WENO scheme 
and the dimension by dimension finite volume WENO scheme can achieve the desired order of accuracy 
while the genuine two dimensional scheme gives smaller errors for the same mesh 

3.3. Oblique shock tubes. The purpose for this test is to see the capability of the rectangular WENO 
schemes in resolving waves that are oblique to the computational meshes For details of the problem, we refer 
to [10] The 2D Sod’s shock tube problem is solved where the initial jump makes an angle 9 against the x 
axis We take our computational domain to be [0, 6] x [0, 1] and the initial jump starting at (x, y) — (2 25, 0) 


9 






Fig 3 1 Oblique Sod’s problem Density p Top contour, genuine two dimensional WENO, middle contour, dimension 
by dimension WENO, bottom cut at the bottom of the computational domain, the solid line is the exact solution, the triangles 
are the genuine two dimensional WENO results, and the circles are the dimension by dimension WENO results 


and making a 6 = j angle with the x axis The solution is computed up to t = 1 2 on a 9G x 1G uniform 
mesh In Fig 3 1 we plot the density contours computed by the above two WENO schemes and the density 
cut at the bottom of the computational domain We can see that both schemes perform equally well in 
resolving the waves The genuine two dimensional scheme gives a slightly better resolution in the contact 
discontinuity and the rarefaction wave 

3.4. A Mach 3 wind tunnel with a step. This model problem is originally from [25] The setup 
of the problem is The wind tunnel is 1 length unit wide and 3 length units long The step is 0 2 length 
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X 


Fig 3 2 Forward step problem, Ax = Ay = jg, jgg, gig from top to bottom 30 contours from 0 12 to 6 41, dimension 
by dimension WENO 

units high and is located 0 6 length units from the left-hand end of the tunnel The problem is initialized 
by a right-going Mach 3 flow Reflective boundary conditions are applied along the wall of the tunnel and 
inflow/outflow boundary conditions are applied at the entrance/exit The corner of the step is a singular 
point and we txeat it the same way as in [25], which is based on the assumption of a nearly steady flow in 
the region near the corner We show the density contours at time t = 4 in Fig 3 2 Only the results from 
the dimension by dimension WENO scheme are shown Uniform meshes of Ax = Ay = i, yig, ^ are 
used 

3.5. Double Mach reflection. This problem is also originally from [25] The computational domain 
for this problem is chosen to be [0, 4] x [0, 1] The reflecting wall lies at the bottom, starting from x = | 
Initially a right-moving Mach 10 shock is positioned at x = g, y = 0 and makes a 60° angle with the x axis 
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Fig 3 3 Double Mach reflection, Ax = Ay — (top and lower left) and (middle and lower right) Genuine two 
dimensional WENO Blow-up regions at the bottom for details 


For the bottom boundary, the exact post-shock condition is imposed for the part from x = 0 to x — ~ and 
a reflective boundary condition is used for the rest At the top boundary, the flow values are set to describe 
the exact motion of a Mach 10 shock We compute the solution up to t — 0 2 Fig 3 3 and Fig 3 4 show 
the equally spaced 30 density contours from 1 5 to 22 7 computed by the genuine two dimensional and the 
dimension by dimension WENO schemes We use uniform meshes with Ax = Ay = 7^, ^ We can see 
that the results from both schemes are comparable 

4. 2D finite volume WENO schemes on triangular meshes. Both third and fourth order finite 
volume WENO schemes on triangular meshes have been constructed in [8] The optional linear weights in 
such schemes are not unique These are then chosen to avoid negative weights whenever possible, and if that 
fails, a grouping (of stencils) technique is used in [8], which works fairly well in the third order case with 
quite general triangulation but can yield positive weights for the fourth order case only with fairly uniform 
triangulation In this section, we do not seek positive linear weights as in [8], but rather use the splitting 
technique to treat the negative linear weights when they appear For scalar equation, the scheme is stable 
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Fig 3 4 Double Mach reflection, Ai = Ay = (top and lower left) and (middle and lower right) Dimension by 
dimension WENO Blow-up regions at the bottom for details 


in all runs For systems of conservation laws, there are still occasional cases of overshoot and instability, the 
reason seems to be related to characteristic decompositions and is still being investigated 

4.1. Accuracy check for a smooth problem. We solve the 2D Burgers equation (1 10) with the 
same initial and boundary conditions as before using the fourth order finite volume WENO scheme [8] The 
solution is computed up to t = when no shock has appeared The meshes used are 1) uniform meshes 
with equilateral tnangulation and 2) random triangulation For the uniform meshes we do not seek positive 
weights as was done in [8], rather we use the splitting technique to treat the negative linear weights when 
they appear Table 4 1 indicates that close to fourth order accuracy can be achieved 

4.2. Discontinuous problem 1: Scalar equation in 2D. Having shown the stable results with the 
splitting treatment of negative linear weights for a fourth order finite volume WENO scheme for the Burgers 
equation in section 2, we now test the fourth order WENO scheme on the Buckley-Leverett problem whose 


13 



Table 4 1 

2D Burgers equation accuracy check 


uniform mesh 

nonuniform mesh 

Ax 

L°° error 

order 

Ax 

L°° error 

order 

2 57E-1 

6 22E-4 


2 67E-1 

2 11E-3 


1 29E-1 

4 61E-5 

3 75 

1 26E-1 

2 35E-4 

2 92 

G 43E-2 

2 18E-6 

4 40 

G 32E-2 

2 90E-5 

3 03 

3 21E-2 

1 38E-7 

3 98 

3 34E-2 

2 G1E-G 

3 78 

1 61E-2 

G 93E-9 

4 32 

1 GGE-2 

2 71E-7 

3 24 

8 08E-3 

6 70E-10 

3 40 

7 44E-3 

1 57E-8 

3 55 


Fig 4 1 2D Buckley- Leverett equation the mesh 



flux is non- convex 


f(u) = 


u 2 + 0 25(1 - u y ’ 


g{u) = 0 


with the initial data u = 1 when —7 < x < 0 and u = 0 elsewhere The solution is computed up to 
t = 0 4 The exact solution is a sliock-rarefaction-contact discontinuity mixture The mesh we use here is a 
non-uniform tnangulation, shown 111 Fig 4 1 Fig 4 2 shows that the waves have been resolved very well 


4.3. Discontinuous problem 2: System of equations in 2D. We consider the 2D Euler equations 
in the domain [—1, 1] x [0, 0 2] The Sod and Lax shock tube initial data is set in the x direction and periodic 
boundary condition is applied in the y direction We use the fourth order finite volume scheme on triangular 
meshes to solve the above problem The mesh we use here is uniform But we do not seek positive weights as 
was done m [ 8 ], rather we use the splitting technique in section 2 to treat the negative linear weights when 
they appear In fact we set deliberately certain linear weights to be negative to test the splitting technique 
Fig 4 3 shows the numerical results of the Sod and Lax problems 

It seems that there are still oscillations and instability for some non-uniform triangular meshes for the 
fourth order WENO schemes applied to Euler equations As the method works well for the same meshes 
with a scalar equation, the problem might be from the characteristic decompositions This is still under 
investigation 


5. Concluding remarks. We have devised and tested a simple splitting technique to treat the negative 
linear weights in WENO schemes This technique involves very little additional CPU time and gives good 
results in most numerical tests The only case where it still yields oscillations and instability is when a fourth 
order finite volume WENO method is used on some non-uniform triangular meshes for Euler equations, the 
reason of which, presumably related to characteristic decompositions, is still under investigation 
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Fig 4 2 2D Buckley-Leverett equation at t=0 f, with splitting Left the solution surface, Right the cut at y = 0 1 (solid, 
line exact solution, symbols numerical solution) 




Fig 4 3 Density plot, Left Sod problem, Right Lax problem, with splitting Roughly 100 points in the x direction 
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